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Abstract 

We propose a new method to calculate electric dipole moments induced by the 
strong QCD 0-term. The method is based on the gradient flow for gauge fields 
and is free from renormalization ambiguities. We test our method by computing 
the nucleon electric dipole moments in pure Yang-Mills theory at several lattice 
spacings, enabling a first-of-its-kind continuum extrapolation. The method is rather 
general and can be applied for any quantity computed in a 0 vacuum. This first 
application of the gradient flow has been successful and demonstrates proof-of- 
principle, thereby providing a novel method to obtain precise results for nucleon 
and light nuclear electric dipole moments. 


28 October 2015 



1 Introduction 


The electric dipole moments (EDMs) of the neutron and proton are very sensitive 
probes of CP-violating sources beyond those contained in the Standard Model (SM). 
In fact, the current bound on the neutron EDM strongly constrains many models 
of beyond-the-SM (BSM) physics. At current experimental accuracies, a nonzero 
nucleon EDM cannot be accounted for by the phase in the quark-mass matrix. 
This implies that such a signal is either caused by a nonzero QCD 6 term or by 
genuine BSM physics which, at low energies, can be parametrized in terms of higher¬ 
dimensional CP-violating quark-gluon operators. Irrespective of the origin, the signal 
for the nucleon EDM will be small and largely masked by strong-interaction physics, 
which presents a formidable challenge to the interpretation of such a signal. To 
disentangle the origin of a nonzero EDM measurement (e.g. 6'-term or BSM), a 
quantitative understanding of the underlying hadronic physics is required. 

The current experimental limit on the neutron EDM is \d]sf \ < 2.9-10“^^ efm [1] and 
experiments are underway to improve this bound by one to two orders of magnitude. 
The bound on the proton EDM is induced from the ^®®Hg EDM limit [2] and is 
\dp\ < 7.9 ■ 10“^^ efm. Plans exist to probe the EDM of the proton directly (and 
other light nuclei) in storage rings [3] with a proposed sensitivity of 10“^® e ■ fm, 
thus improving the current bounds by several orders of magnitudes and covering a 
wide range where BSM physics can show its footprint. 

Nucleon EDMs arising from the QCD 0-term or BSM physics have been calculated 
both in models [4] and in chiral perturbation theory [5,6]. In the latter approach, the 
nucleon EDMs are calculated in terms of effective CP-odd hadronic interactions that 
have the same symmetry properties as the underlying CP-odd sources at the quark 
level (for a review, see [7]). The calculated EDMs depend on several low-energy 
constants (LECs) whose sizes are in most cases unknown and need to be estimated 
or calculated with lattice QCD. 

Lattice QCD can thus be used to perform an ab initio calculation of the nucleon 
EDM. For the 6*-term, this has already been shown in the pioneering works in 
refs. [8,9] and later in ref. [10] (for BSM sources only the nucleon EDMs arising 
from the quark EDMs have been calculated with lattice QCD [11]). The chiral and 
inhnite volume extrapolations of unpublished lattice data from Shintani et ah have 
been performed in refs. [12,13]. The calculation of the EDM within a lattice (dis¬ 
cretized) formulation of QCD is very non-trivial, and present large difficulties for 
two main reasons. The renormalization of the CP-odd operators and the degradation 
of the signal-to-noise ratio towards the chiral limit. Additionally, the 6'-term itself 
introduces an imaginary term in the real Euclidean action, which produces a sign 
problem and precludes the use of standard stochastic methods employed by lattice 
QCD. Ref. [14] performed a lattice QCD calculation of the neutron EDM induced 
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by a 0-term that was analytically continued into the complex plane. This allows the 
usage of standard stochastic methods. 

In this paper we propose, without relying on any complex rotation of the 6'-term, a 
method based on the gradient flow for the gauge fields [15] that has no renormaliza¬ 
tion ambiguities and, to our knowledge, is the only method that allows a theoretical 
sound continuum limit. A hrst account of this method can be found in ref. [16]. 

The remainder of the paper is organized as follows: The next section gives a cur¬ 
sory discussion of the phenomenology of the nucleon EDM. In sec. 3 we introduce 
dehnitions and our method. Sec. 4 discusses the gradient flow for gauge helds and 
its relevance to the calculations presented in this paper. We provide details of our 
lattice calculations and their results in secs. 5 and 6, followed by a discussion in the 
ensuing section. 


2 Phenomenology of the QCD theta term. 


The discrete space-time symmetries parity P, time-reversal T, and the combination 
of charge conjugation and parity CP, are broken in QCD by the QCD 6 term. In 
the case of three quark flavors the QCD action is given by 

Se = J d'^x [£qcd - iOq{x)] , (2.1) 

where Dqcd is the standard Euclidean QCD Lagrangian 

£qod = + M)ii ( 2 . 2 ) 

and 

Qi^) = ^:^^t^upaF^Ax)FpAx ), (2.3) 

is the topological charge density. The fermion held containing up, down and strange 
quarks is denoted by A = {u ,d, s)^ and is the gluon held strength tensor. 

(eoi 23 = +1) is the completely antisymmetric tensor, the gauge-covariant 
derivative, M the real 3x3 quark-mass matrix, and 6 the coupling of the CP- 
odd interaction. In eq. (2.2) the complex phase of the quark-mass matrix has been 
absorbed in the physical parameter 6, i.e. we choose a fermionic basis where the 
CP-odd interaction comes solely from the topological charge density. 

The most important consequence of the QCD 6*-term is that it induces EDMs of 
hadrons and nuclei. The hrst dedicated EDM experiment was the neutron EDM 
experiment in 1957 [17]. Since then, the accuracy of the measurement has been 
improved by six orders of magnitude without hnding a signal. The current bound 
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(In < 2.9-10 fm [1] sets strong limits on the size of 6 and sources of CP violation 
from physics beyond the SM. 


In order to set a bound on the 6 term, it is necessary to calculate the dependence of 
the neutron EDM on 6. One way to do this is by using yPT. This calculation has 
been done up to next-to-leading order (NLO) in both SU{2) [18,19,6] and SU{3) 
[21,20,5] yPT. Focusing here on the two-flavored theory, the neutron (d^v) and proton 
EDM [dp) are given by: 


dN — 


egA9o 

167r2F2 


In 


M2 




^Af,EDM 


2M, 


N, 


(2.4) 


dp — 




In 


Ml 


2'kMt, 


(2.5) 


167r2F2 Ap,EDM , 

in terms of gA — 1-27 the strong pion-nucleon coupling constant, i^r — 92.2 MeV 
the pion decay constant, and Mn the pion and nucleon mass respectively, e > 0 
the proton charge, and, in principle, three low-energy constants (LECs) of CP-odd 
chiral interactions: g^ and (iN,p- The first one, g^, is not free as discussed below 
(see eq. (2.9)). The latter two are absorbed in renormalization-scale, /i, independent 
constants 

STT^F^dNilA) 


Aat.edm — exp 


e^Ado 


f 8Tr^F^dp{g,) 

Ap,EDM = H exp < - — - 

I egAgo 


( 2 . 6 ) 


(2.7) 


The hrst term in brackets in eqs. (2.4) and (2.5) arises from the leading-order one- 
loop diagram involving the CP-odd vertex 


^■KN{d) = 


Ntt-tN, 


2P 


( 2 . 8 ) 


in terms of the nucleon doublet N and the pion triplet tt. The LO loop is divergent 
and the divergence and associated scale dependence have been absorbed into the 
counter terms dv.p which signify contributions to the nucleon EDMs from short- 
range dynamics and appear at the same order as the LO loop diagrams. The second 
term in brackets in eqs. (2.4) and (2.5) is a next-to-leading-order correction. 


The 6*-term breaks chiral symmetry as a complex quark mass. As such, chiral sym¬ 
metry relates to known CP-even LECs [21,22]. In particular, it is possible to write 
[23] 


where {Mn — 
splitting, rh = 


go 


(M, 


N 


Mp) m* 


Me = (15.5 ±2.5) ■ 10 
me 


-3 


9 


(2.9) 


2P^ 2P^ 

Mp)s*™ng jg quark-mass induced part of the proton-neutron mass 
{mu ± md)/2,mi, = mumdj{mu ± md ), and e = (m„ - md) /(m„ ± md) ■ 
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To get a rough estimate of the sizes of the nucleon EDMs, we can insert eq. (2.9) in 
eqs. (2.5) and (2.4). If we assume that Aedm — ATat, we obtain 

~-2.1 ■ 10-^0 efm, (2.10) 

dp ~+2.5-lO“^0efm, (2.11) 

as a rough estimate of the nucleon EDMs. A comparison with the experimental 
bound then gives the strong constraint 6 < 10“^°. Clearly, a more reliable constraint 
on 6 requires a hrst-principle calculation of the nucleon EDMs. In the isospin limit, 
Qq scales as m ~ such that the loop contributions to the EDMs vanish in the 
chiral limit as M^logM^ (see eqs. 2.4 and 2.5). 

In the isoscalar combination + dp the loop contribution cancels out to a large 
extent. For observables sensitive to this combination, such as the deuteron EDM 
[24,25], a hrst-principle calculation of the total nucleon EDM is important to differ¬ 
entiate the d-term from possible BSM sources of CP violation [26,27]. In the specihc 
case of the isoscalar combination a precise evaluation of disconnected diagrams is 
needed in any lattice QCD calculation. 


3 The electric dipole moment 


The theory is dehned in Euclidean space with the action given in eq. (2.1). The 
EDM of a nucleon is related to the spatial charge density distribution. If we dehne 
the quark charges as Qu = 2/3 e and Qd = Qs = —1/3 e, the nucleon EDMs are 
obtained from the matrix element of the electromagnetic current 

J^{x) = Quu{x)')^u{x) + Qdd{x)')^d{x) -h Qss{x)'^^s{x ), (3.1) 

between nucleon states in the 9 vacuum, 

{N\v\s')\J,\N\v,s))=u%{v\s')T,{Q^)u%{v,s). (3.2) 

has the most general four-vector structure consistent with the symmetries of 
the action (2.1) such as gauge, 0(4), C and CPT invariance . Note that the photon 
momentum q = p' — p in Euclidean space is 

Q, = (^4, Q) = q) , = -{q^ + |qP = -q^ • (3.3) 

Following ref. [8] the dependence of the matrix element is parametrized by a linear 
combination of CP-even and CP-odd form factors. Using Euclidean 0(4) rotational 
invariance, gauge symmetry and the spurionic symmetry P x 0 —)■ —6, the most 
general decomposition reads 

r,(Q^) = g{9^)r<;:^\Q^) + teh{9^)Tf\Q ^), (3.4) 
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where g, h are even functions of 9. The CP-even contribution is given by 

rr”(Q') = 7 ^-Fi(Q") + (3.6) 

where the Dirac and Pauli form factors Fi and F 2 are related to the electric and the 
magnetic form factors 

Ge(Q^)=F,(Q")-^F 2(Q"), G„(Q2)=f,(Q2) + f2(Q"). (3.6) 

The CP-odd term reads 

T‘‘(Q") = ff;,.75^F3(Q") . (3.7) 

In the literature P°‘^'^(Q^) usually contains an additional parity violating form factor, 
the anapole form factor. The anapole form factor breaks parity symmetry but does 
not break time reversal, i.e. is CP-even while breaking both C and P. It therefore 
does not contribute to the amplitude in eq. (3.2) of the electromagnetic current 
evaluated in a 0 vacuum. In other words, the ^-term alone cannot induce a nucleon 
anapole form factor. The EDM is directly related to the CP-odd F^{Q‘^) form factor 
at zero momentum transfer 


dN = 0g{9^ 


Fi^iO) .F 


2M, 


~ 9 


3^(0) , ^ 7 . 3 , 


N 


2M. 


+ o(r 


N 


(3.8) 


In lattice calculations, matrix elements can be extracted from the large-distance 
behavior of appropriate correlation functions in Euclidean space-time. In the case 
at hand, one considers three-point correlations such as 

~ , (3-9) 

where Af is an interpolating operator with the same quantum number of the nucleon. 
The three-point functions are to be evaluated with the Euclidean action Sq. Although 
the action with 9^0 cannot be directly studied by numerical Monte Carlo methods, 
in the small 9 limit one can obtain the desired result for EDM by expanding around 
9 = 0 and taking only the linear term in 9. That is, for a generic expectation value 
of product of operators, O, in a 6*-vacuum, we can write 

{0)e ~ {O)e=o + ^0{O J d^x q{x))e=^ + 0 ( 0 ^), ( 3 .IO) 

where q{x) is the topological charge density (2.3). The experimental bound on 9 is 
currently 9 < 0(10^°) (see sec. 2), thus a power expansion in 9 is well justihecQ 

^ Alternatively, the nucleon EDM at finite 6 can be also determined using reweighting 
techniques with the complex weight factor 
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Fig. 1. Flow-time dependence of the topological charge for two different gauge fields. 

In general this proposal could be hampered by the impossibility of giving a sound or 
practical dehnition on the lattice of the topological charge density and its continuum 
limit. In this work we propose to directly compute the matrix element 

{O J q{x))g=o , (3.11) 

using the gradient flow (see sect. 4) to dehne the topological charge density [15]. By 
doing so, we have a theoretically sound dehnition of the correlation function with 
no renormalization ambiguities and a well-dehned continuum limit. 


4 Gradient flow and the topological snsceptibility 


The gradient flow [15] of Yang-Mills gauge fields is defined as follows 

dtBf^ = 1 (4-1) 

where the flow-time t has a time-squared dimension, 

Gf_iu = — dyBfj^ + [5^, By ], + [5/^, ■ ], (4-2) 

and the initial condition on the flow-time-dependent field B^{t, x) at t = 0 is given by 
the fundamental gauge field. The gradient flow for gauge fields and for fermions [28] 
has several applications and here we mention the definition of a relative scale [15,29], 
the determination of the strong coupling constant [15,30] and of the chiral conden- 
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sate [28,31], the calculation of the energy-momentum tensor [32,33] and of the topo¬ 
logical susceptibility [34,35]. We have recently proposed to use the gradient flow for 
the determination of the strange content of the nucleon [16]. 

One way to understand the flow equations is to consider them as steepest descent 
equations in the space of gauge helds. As such the evolution along the flow drives 
the gauge conhgurations towards local minima of the action. The topological charge 
is dehned at non-vanishing flow-time as 

Q{t) = J d^x q{x,t) , (4.3) 

with 

Qix, t) = t). (4.4) 

In £g. 1 we show the flow-time evolution of Q{t) evaluated on two representatives 
of our gauge ensembles. Rather rapidly Q{t) reaches a plateau where it assumes an 
almost integer value saturating the corresponding instanton bound. 

Another way to understand the effect of the gradient flow on the gauge fields is 
apparent already at tree-level. The smoothing at short distances over a range ^/^ 
corresponds in momentum space to a Gaussian damping of the large momenta. This 
results in a very interesting property of the flowed gauge fields B^{x, t): they are free 
from ultraviolet divergences [15,36] for all f > 0 and do not require any renormal¬ 
ization. This powerful result can be used to simplify the renormalization pattern of 
operators involving gauge fields. In general, one would need to relate the local opera¬ 
tors evaluated at non-vanishing flow-times with the ones at zero flow-times. The case 
of the topological charge and all the correlation functions containing the topologi¬ 
cal charge is special, because in this case we can define the topological charge, and 
for example the topological susceptibility, directly at non-vanishing flow-time [15]. 
The Euclidean theory is prepared on a 
lattice of spacing a and volume X T. 

The calculations in this and the follow¬ 
ing sections have been performed with 
the standard Wilson gauge action, with 
f3 = at 4 different lattice spacings 
a = 0.093,0.079,0.068,0.048 fm corre¬ 
sponding to /3 = 6.0, 6.1, 6.2, 6.45. In this 
work we use the Sommer scale [37,38], 
ro = 0.5 fm, to fix the lattice spacing 
in physical units. The size of the box 
is respectively L/a = 16,24,24,32 with 
T/L = 2. To generate the gauge links, we have used a heat bath algorithm with a 
ratio of number of over-relaxation steps, Nqr-, over a number of Cabibbo-Marinari 
updates. Nub, per sweep of Nqr/Nhb = 4/1. For thermalization we have per- 


/3 

iVth 

iVup 

N, 

N 

’ meas 

6.0 

2000 

200000 

1000 

1000 

6.1 

2000 

65000 

325 

325 

6.2 

2000 

60000 

300 

300 

6.45 

2000 

122400 

612 

153 


Table 1. Summary of our runs: is the 

number of thermalization updates; N^p 
is the total number of updates; Ng is the 
number of gauges saved and A^meas is the 
number of gauges analyzed. 






Fig. 2. Left plot: flow-time dependence of the topological susceptibility at several lattice 
spacings. Right plot: continuum limit of the topological susceptibility. The yellow band is 
a linear extrapolation in a? compared with a constant fit. 

formed 2000 updates. For the finest lattice spacing, we have analyzed all correlation 
fnnctions skipping 800 ganges while for the remaining correlation fnnctions we have 
skipped 200 ganges. A snmmary of parameter rnns is given in tab. 1. 


With these choices we have observed no signihcant antocorrelation for all onr lattice 
spacings. We obtain the same ontcome also for the correlation fnnctions nsed for 
the determination of the EDMs. A more detailed discnssion of antocorrelations for 
the fermionic correlation fnnctions is given in sec. 5. 


The gradient flow eqnation at hnite lattice spacing is solved following app. C of 
ref. [15] with step-size for the flow-time e = 0.01. The topological charge density 
is dehned as in eq. (4.4) where is the lattice implementation of the held 

tensor dehned in ref. [39]. Any other dehnition of the topological charge density in 
a pnre Yang-Mills theory reqnires a hnite mnltiplicative renormalization [40] that 
has to be determined as a fnnction of the bare conpling, in order to perform the 
continnnm limit. With the dehnition based on the gradient how, this renormalization 
factor is 1 independently of the lattice action nsed. Using the dehnition of the 
topological charge density at non-vanishing how-time given above, we can also dehne 
the topological snsceptibility 

Xt{t) = ^ j d^x d^y{q{t,x)q{t,y)) . (4.5) 

The topological snsceptibility dehned as in eq. (4.5), bnt at vanishing how-time 
t = 0, not only needs a mnltiplicatively renormalization, bnt, more importantly, has 
a 1/a^ power divergence However, with the dehnition at non-vanishing how-time, 

^ A notable exception is the definition of the topological susceptibility proposed in ref. [41] 
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the topological susceptibility needs no renormalization and it has a well dehned 
continuum limit. We have computed Xt for several lattice spacings as a function 
of the flow-time. In the left plot of fig. 2 we show the topological susceptibility in 
physical units as a function of v^/ro- The divergence-free property of the gauge 
helds at non-vanishing flow time allows us to perform the continuum limit at hxed 
value of v^/ro- In the continuum limit we expect the topological susceptibility to 
be flow-time independent for every positive flow-time, v^/^o > 0 [15]. 

For small flow-time values we observe two different effects. First, for < 0.1 
fm we observe a rapid increase of Xt that is just a reflection of the short distance 
singularities discussed above. Second, for 0.1 fm < < 0.2 fm we observe some 

discretization effects. For y^ > 0.2 fm we hnd complete agreement between all 
lattice spacings and, as expected, Xt is flow-time independent. We perform the con¬ 
tinuum limit at y^/ro = 0.8 and this is shown in the right plot of £g. 2 where we 
compare a linear extrapolation in a? with a constant one. We decide to quote as 
final result 

1/4 = 195.9(4.9) MeV, (4.6) 

that is a constant fit including all lattice spacings. The values at all lattice spacings 
and different extrapolations to the continuum limit are listed in tab. 2. This result 
is in perfect agreement with the result [42] obtained using the index theorem with a 
chiral lattice Dirac operator and the result [43] obtained using the spectral projector 
method. Very recently a paper has been submitted with a precise determination of 
the topological susceptibility [44] using the gradient flow. The results are consistent 
within statistical uncertainties. 


5 CP-broken vacuum and nucleon mixing 


The form factor F 3 , directly related to the nucleon EDM, defined in eq. (3.7) can 
be computed non-perturbatively with suitable ratios of the following 2- and 3-point 
functions in a 6 * vacuum 

<T?,r^(p,a;o) = a^[^e*P^(A/'(x,a;o)7^(0))^ , (5.1) 

X 

G‘„j^„(VuV 2 ,x„,y„) = (A/'(x,i„)J,(y,j/„)V(0))^ . (5.2) 

x,y 

Here, the baryon interpolating fields are 


J\f{x) = SABCUAix) 


u%{x)C'y5dc{x) , 


(5.3) 


based on spectral projectors. 
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JI{x) = Eabc u^{x)C'y 5 (fB{x) Uc{x), (5.4) 

and C is the charge conjugation matrix. We now describe in some detail the spectral 
decomposition for the 2-point functions and defer to app. A for the slightly more 
cumbersome spectral decomposition of the 3-point functions. Most of the discussion 
on the spectral decomposition follows Shintani et al. [8], but we rederive some of their 
results for clarity and to be consistent with our normalizations. The key ingredient 
of the spectral decompositions is the matrix element of the interpolating operator 
of the nucleon between the 6 vacuum and a single nucleon state 

(0|A/'|A^®(p,s)) = Zn{0)u%{p,s) . (5.5) 

In a theory that does not preserve parity, for instance due to the presence of a 6- 
term, the nucleon state does not have a dehnite parity and the nucleon spinor can 
be written as 

u%{p, s) = UNip, s) , (5.6) 

where U]\f{p,s) is the nucleon spinor in the 6 = 0 vacuum. In other words, the 
nucleon spinor satishes the modihed Dirac equation 

+ Mn{9) = 0 . (5.7) 

The theory still preserves the spurionic symmetry Pq : P x 6 ^ —6, where P is 
the standard parity transformation. This implies that both the energies and the 
amplitudes M{6), Z{6) are even functions of 9, M{9) = M + 0{9‘^) and Zj^ifi) = 
Zn + 0{9^). 

The phase a]sf{9) plays a very important role in the determination of the EDM. From 
the spurionic symmetry Pq we deduce that aN{9) = —aN{—9) and for small values of 
9, aN{9) = a^^^9+0{9^). It is important to determine precisely the mixing parameter 
ajsf before extracting the CP-odd form factors from the 3-point functions. The reason 
is that the mixing between different parity states can induce a spurious CP-odd 
contribution to the correlation function proportional to the CP-even form factors. 
These spurious contributions can be subtracted only with a precise determination 
of the mixing angle aN{9). The details of these spurious contributions and relative 
subtractions are detailed in app. A. 

For on-shell nucleons with energy —ipo = E]\f{p) where E]sf{p) = \J\p\^ + M^, the 
inhnite volume normalization reads 

{N\q,s)\N\k,s')) = (27r) V2F^(0; q)^2E^{9; k)<5(3)(k - . (5.8) 

Taking into account the parity mixing, the completeness relation of the nucleon 
spinors with spatial momentum p reads 

u%{p, s)u%{p, s) = En{9-, p) 7 o - i'jkPk + . (5.9) 
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For small values of we have 

J2u%ip,s)u%{p,s) = ^jv(p)7o - ilkPk + Mn (l + + 0(9^). (5.10) 

5 

We can now perform the spectral decomposition of the nucleon 2-point functions in 
a 0-vacuum. Retaining only the one-state leading contribution we obtain 

-EN{e-,p)xo 

G%n{p,xo)= , \Zn{0]p)\ ^m^(p,s)m^(p,s), (5.11) 

ZhN[9;p) g 

and using the completeness relation we get 

p-En{S-,p)xo 

G%N{P,xo)a0 = in —7l^^(^;p)l ^iv(0;p)7o-mPfc + Mjv(0)e2*“^(^)^5 

Zt/N\P]P) -'"P 

(5.12) 

where a and (5 are the Dirac indices. Expanding the l.h.s. of eq. (5.12) in powers of 
0, we obtain 

G%j^{p,Xo) = Gmn{p, Xq) +ieG%pf{p,Xo)+ 0{e‘^) (5.13) 

where 

G'Ar7v(p,a;o) = ^ e*P^ ^A/'(x, a;o)W'(0)) , (5.14) 

X 

and 

Gnn{p,xo) = (N'{-k,xo)W{())Q) . (5.15) 

X 

The term linear in 0 can be computed inserting the topological charge in the nucleon 
2-point function. The topological charge Q, defined in eq. (4.3), is computed as 
detailed in sect. 4 using the gradient flow. In this way the topological charge is 
free from any renormalization ambiguity and the continunm limit can be safely 
performed keeping hxed the flow-time in physical nnits. To minimize discretization 
effects we choose, \/^/ro = 0.8. We omit the flow-time dependence of Q{t) because 
in this range of flow-times any correlator involving the topological charge is flow 
time independent (see sect 4). 

By expanding the spectral decomposition, i.e. the r.h.s. of eq. (5.12), in powers of 
0, we obtain the standard nucleon spectral decomposition 

g-£^]v(p)a;o 

GAr7v(p,a;o) = , . I^w(p)| [^v(p)7o - ilkPk + M] , (5.16) 

ZrjN\P) 

and the term linear in 0 

p-En{p)xo 

G%n{p^xo) = ^ |2iv(p)| 2Mjvaff75- (5.17) 

ZXjn{P> 
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Fig. 3. Distribution of the topological charge computed for \/^/ro = 0.8 for 4 different 
lattice spacings. 

For simplicity we have not written down the opposite parity states propagating from 
T. If we project to p = 0 and to positive parity we obtain 

C{xo) = tT[P+GMO,Xo)] = 2\Z^\^e-^-^° + --- , (5.18) 


and 


C^{xo) = tr [P+75G^^(0,a;o) 


= 2 




(5.19) 


We observe that the two correlators have the same leading exponential behavior. If 
the sampling of all topological sectors is correctly performed, the effective masses of 
the two correlators should agree asymptotically for large Euclidean times. 


In Eg. 3 we plot the distribution of the topological charge for 4 different lattice 
spacings at = O.SrQ. Details on the definition can be found in sec. 4. The 
distribution looks reasonably Gaussian with all average values statistically consistent 
with zero. We observe for a = 0.079 fm, that the distribution has slightly larger 
width, but this is related to the slightly larger physical volume of that lattice. As we 
have seen in the previous section, the topological susceptibility does not show any 
sign of discretization errors. 
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■ 10^ [fm^ 


x]'^ [MeV] 

OtN 

8.675135 

198.1(2.3) 

0.289(19) 

6.250475 

199.7(4.8) 

0.314(38) 

4.615747 

195.1(4.3) 

0.324(33) 

2.285814 

191.7(5.9) 

0.301(35) 

0 [fit 1] 

195.9(4.9) 

0.314(35) 

0 [fit 2] 

191.4(6.0) 

0.326(40) 


Table 2 

Numerical results for the topological susceptibility and the CP-odd mixing angle for 
several lattice spacings. The continuum extrapolated values are obtained with a constant 
fit using the 3 finest lattice spacings (fit 1) and a linear extrapolation in using all the 
lattice spacings (fit 2). 

For the computation of the 2-point functions, we have studied 3 different levels 
Si, i = 1,2,3 of Gaussian smearing [45]. The relevant parameters of the Gaussian 
smearing that we have considered, usually labeled as {a,NG}, are si = {2,30}, 
S 2 = {4,25}, S 3 = {5.5,70}. We have found that the S 3 smearing has a better 
projection on the fundamental state, but it is also the smearing that adds more 
noise to the correlator. Gompromising between an earlier plateau and a less noisier 
correlator, we have decided to choose the smearing S 2 for the 2 coarsest lattice 
spacings and the smearing S 3 for the 2 hnest spacings. 

The fermion lattice action is the non-perturbatively improved Wilson action [46,47,48] 
The propagators are computed with sources located stochastically in the 3 spatial 
directions. We choose 20 stochastic spatial points for the hnest lattice spacing and 
10 stochastic points for the others. The rational behind this choice is to have 0{L/a) 
different stochastic points to improve the overlap between the topological charge and 
the fermionic part of the correlation functions. We stress that this is very important 
to improve the signal-to-noise ratio not only of the 2 -point functions, but especially 
for the 3-point functions which we discuss in the next section. 

We have performed the calculation at 4 lattice spacings (see sec. 4) and at the 
following set of momenta 

27r 

{P} = — ■ {(0, 0, 0), (±1, 0, 0), (±1, ±1, 0), (±1, ±1, ±1), (± 2 , 0, 0)} . (5.20) 

The values of the quark mass for all the lattice spacings corresponds to a value of 
the pseudoscalar mass, Mps — 800 MeV, hxed in physical units [49]. From coarser 
to hner spacings they correspond to the following values of the hopping parameter, 
K = {0.13353,0.13423,0.13460,0.13485}. For these values of n we have computed 
the nucleon mass, that shows very small discretization errors and it corresponds to 
a value ~ 1.65 GeV. We have also checked the dispersion relation for all lattice 
spacings and hnd that the discretization errors are below our statistical accuracy. 
In £g. 4 we show our results for our coarsest lattice spacing of the nucleon energy 
squared for all the |pp of the set (5.20) with the continuum form of the dispersion 
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Fig. 4. Results for the nucleon energy squared at different values of |pp compared with 
the continuum dispersion relation E‘j^ = + |pp. 

relation. 


In the left plot of fig. 5 we show the effective masses of the two correlators in 
eqs. (5.18,5.19) for the finest lattice spacing. It is clear that for large Euclidean 
times we have perfect agreement between the two effective masses and very similar 
results are obtained for all the other lattice spacings we have. This is also conhrmed 
on the right plot of fig. 5 where we show, again for our finest spacing, the nucleon 
mass obtained from the 2 correlators in eqs. (5.18,5.19) for different fit ranges. The 
calculation of the mixing angle is now straightforward 


tr 


( 0 , ^o) 


tr [P+Gjyjv ( 0 , Xo)] 


^ H- 


(5.21) 


and we expect a plateau for large Euclidean times with higher-order corrections 
that are exponentially suppressed. In the left plot of fig. 6 we show the Euclidean 
time dependence of obtained from the ratio in eq. (5.21) at a = 0.048 fm. A 
plateau is easily identified as is the case for all the other lattice spacings. This is 
just a reflection of the previous result, namely that asymptotically both correlators 
in eqs. (5.18,5.19) are dominated by the same exponential behavior with the same 
mass. 


We have performed several checks on the calculation of the mixing angle because 
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Fig. 5. Left plot: comparison of the effective masses in lattice units obtained from the 
nucleon correlators at a = 0.048 fm with and without the insertion of the topological 
charge, i.e. from eq. (5.18) (Mgff) and from eq. (5.19) (M^). Right plot: comparison of 
the nucleon masses in lattice units obtained with the nucleon correlators (5.18), Mat, 
and (5.19), M®, at o = 0.048 fm for different fit ranges = 28). 

a solid determination of is crucial for a correct and precise extraction of the 
nucleon EDM as detailed in the next section and in the app. A. We can determine 
the mixing angle from ratios as in eq. (5.21) but with correlators projected at non¬ 
vanishing spatial momenta. If we choose the same interpolating operators for the 
two correlators, from the spectral decomposition in eqs. (5.14,5.15) we obtain 


Ejv(p) RMjy tr [P+'y5G%NiP,Xo) 
2Mm tr [P+G'wAr(p, xo)] 


_ ^,( 1 ) 


= a 


N 




(5.22) 


Up to discretization effects, should not depend on the momentum chosen in the 
nucleon 2-point functions. In the right plot of £g. 6 we show the |pp dependence of 
for our hnest lattice spacing. We expect an increased uncertainty as we increase 
the nucleon momentum and we see perfect agreement between all the values of the 
mixing angle. We obtain similar results for all the other lattice spacings. 


Another check of our calculation concerns the autocorrelation time of the correlators 
containing the topological charge. Critical slowing down has been observed for topo¬ 
logical charge and susceptibility both in QCD and Yang-Mills theory [51] using an 
Hybrid Montecarlo (HMC) algorithm. In particular, it is expected that the problem 
can become relevant for lattice spacings below 0.05 fm. Even though the situation 
here is different because the correlators contain explicitly fermionic propagators and 
we do not use an HMC algorithm, at our hnest lattice spacing a = 0.048 fm we 
have computed the autocorrelation function for all the Euclidean times xq of the 
correlator G®^(0,xo) and We have followed refs. [50,52] for the determination 
of the autocorrelation function and integrated autocorrelation times. 
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Fig. 6. Left plot: Euclidean time dependence of determined from the ratio in eq. (5.21) 
at o = 0.048 fm. The red band indicates our choice for the plateau and the corresponding 
uncertainty. Right plot: momentum dependence of as determined from the ratio in 
eq. (5.22), at a = 0.048 fm. 


normalized autocorrelation of a{Xp/a=10) 



x„/a 




Fig. 7. Left plot: normalized autocorrelation function and estimate of the integrated au¬ 
tocorrelation time Tint with the automatic windowing procedure [50] of at xq/u = 10. 
Right plot: estimate of Tint of for all Euclidean times. The dashed red line indicate 
the absence of autocorrelation, Tint = 0.5. 

In the left plot of fig. 7 we show for aX Xq/ a = 10 the normalized autocorrelation 
function and the estimate of the integrated autocorrelation time Tint with an auto¬ 
matic windowing procedure [50]. On the right plot we show the estimate of Tint of 
for all Euclidean times. It is clear there is almost no autocorrelation. Tint = 0.5, 
for all Euclidean times. As a further check, we compared the error estimate of 
using the autocorrelation function method to a standard bootstrap method. This is 
shown in Eg. 8 from which it becomes clear that we can safely use a boostrap analysis 
to determine our statistical uncertainty for all correlators at all our lattice spacings. 


17 













































































0.4 


0.3 


0.2 


a;;'(x„/a) Bootstrap 


T. 


int 


in 


Xo/a 


Fig. 8. Comparison of error estimates using the autocorrelation function method, labeled 
by Tint, and a standard bootstrap method. 


We can now perform the continnnm limit of for fixed valne of the pion mass. In 
fig. 9 we show the continnnm limit and in tab. 2 we list all the valnes at all lattice 
spacings and different extrapolations to the continnnm limit. We compare a linear 
extrapolation in (yellow band) with a constant extrapolation inclnding the three 
hnest lattice spacings. We observe a perfect agreement for all the extrapolations and 
tiny discretization errors. The theory is non-pertnrbatively improved so we expect 
an O(a^) scaling behavior. 


Since we see no signs of discretization errors, as a hnal resnlt we qnote the valne 
obtained nsing a constant £t exclnding the coarsest lattice spacing 


= 0.314(35). 


(5.23) 


We stress that this is the hrst time that a continnnm limit is performed for this CP- 
mixing angle. The normalization chosen for and the convention for the Dirac 
7 matrices is consistent with the one of ref. [9]. Onr resnlt in the continnnm limit 
differs by 2a from the resnlt of ref. [9] that is obtained with a different fermionic 
and gange action, at a single lattice spacing of a ~ 0.15 fm and at a similar qnark 
mass valne. 
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Fig. 9. Continuum limit of The yellow band is a linear extrapolation in and it is 
compared with a constant extrapolation including all the lattice spacings and excluding 
the coarsest one. 


6 Nucleon electric dipole moment 


The spectral decomposition of the 3-point functions 

G%j^NiPuP2,xo,yo) = (A^(x,a:o)J^(y,l/o)A/'(0))^ , (6.1) 

x,y 

relevant for the determination of the nucleon EDM is detailed in app. A and the hnal 
result for the leading exponentials is given in eq. (A.4). By taking suitable ratios of 
2- and 3-point functions we can extract CP-even and CP-odd form factors dehned 
in eqs. (3.5,3.7). 

The 3-point functions have been computed for the set of momenta {P} dehned 
in eq. (5.20) and, when possible, we have averaged over all equivalent momenta 
conhgurations. We have tested several sink locations xq and, after some numerical 
experiments, we have chosen the following set 

{xo} = (16a, 20a, 20a, 28a), (6.2) 

from the coarsest to the hnest lattice spacing. We work in the SU(3) havor-symmetric 
limit such that the disconnected contributions vanish. 
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Fig. 10. Euclidean time dependence of the proton electric form factor determined from the 
ratio (A.21) for 2 different values of Q^. 


For all the form factors calculations we use a local vector current. The normalization 
constant Zvig^) is taken from ref. [53]. To compute the EDM we use the ratio of 
eq. (A.27). In order to determine we need to subtract contributions propor¬ 

tional to the mixing angle and the CP-even form factors GEiQ"^) and Gm(< 5^)- 
To determine and Gm^Q"^) we use the ratios in eqs. (A.21) and (A.26). In 

fig. 10 we show 2 typical plateaus for the electric form factor of the proton for 2 
different momenta. A plateau is easily identified as for all the other momenta and 
lattice spacings. We collect in a single plot, fig. 11, the electric form factor Ge(Q^) 
of the proton and the neutron including all the lattice spacings. Discretization errors 
are well under control. The red curve for the proton electric form factor is a fit to 
the lowest 4 momenta using the standard dipole form 

^ (6.3) 

b + Sli 


20 

















Fig. 11. Momentum dependence of the electric form factor for the proton (left plot) and 
the neutron (right plot). The red curve for the proton form factor is a phenomenological 
fit restricted to the 4 lowest momenta. The fit result is = 1.5182(71) GeV^ (see text). 

The neutron electric form factor vanishes at = 0 and is rather small for larger 
values of Q^, but we are still able to identify a clear signal over a wide range of 
Q^. This will be important for the determination of the neutron EDM. In Eg. 12 we 
show the same result for the magnetic form factors and the same type of dipole £t 




f^P,N 


1 + 


Q2 


2 ’ 


(6.4) 


where the anomalous magnetic moments kp^n are hxed to their phenomenological 
values. The £t parameters M|, are not expected to reproduce the phenomenolog¬ 
ical values (see for example ref. [54]). We perform these hts as a check that our 
lattice data can be htted by a dipole form, however we do not use these hts in the 
determination of the EDMs below. 

With the precise determination of the CP-even form factors and the mixing angle 
we can now determine the nucleon EDM. By evaluating the ratio on the l.h.s. 
of eq. (A.27) and subtracting the spurious contributions from the r.h.s, we can 
determine the CP-odd form factor F^{Q‘^). In £g. 13 we show the plateau for the 
normalized CP-odd form factor F3(Q^)/2Mjv for all the momenta at our hnest lattice 
spacing. For the lowest 2 momenta it is possible to extract a signal while for the 
largest 2, the signal is too small and consistent with zero. In £g. 14 we show the 
dependence of /2M]si on a single plot including all our lattice spacings. The 

different lattice spacings results agree rather well, indicating small discretization 
errors within our statistical accuracy. For this reason, we determine the EDM by 
extrapolating using the three finest lattice spacings result to = 0 with the same fit 
function. An extrapolation using all four lattice spacings give completely consistent 
results as shown in Eg. 14. For the extrapolation to = 0, we use the SU(2) xPT 
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Fig. 12. Same as fig. 11 but for the magnetic form factors. The fit results are 
= 0.5621(20) GeV^ for the proton and = 0.5224(15) GeV^ for the neutron. 
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Fig. 13. Euclidean time dependence of F3(Q^)/2Mjv in e • fm extracted from the ratio in 
eq. (A.27) for the 4 non-vanishing lattice momenta at our finest lattice spacing. 
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Fig. 14. Momentum dependence of the CP-odd form factor /2Mjq of the proton 

and the neutron. The yellow band is a linear extrapolation in as suggested from yPT 
at NLO using all four lattice spacings. Different colors represent different lattice spacings: 
(5 = 6.0 (cyan), (3 = 6.1 (red), j3 = 6.2 (magenta), j3 = 6.45 (green). As a comparison we 
plot, slightly displaced, the = 0 extrapolation using the three finest lattice spacings. 

result of [6] as a guideline. There the form factor is expanded as 

pP/N/ q2\ 

——-= dp/isf + Sp/nQ'^ + Hp/NiQ"^) ■ (6.5) 

ZIVlpi 

The values at = 0 are the nucleons EDMs and the slope in at small 
Sp/N, are the so called Schiff moments [55]. The functions Hp/n{Q^), defined in [6], 
scale as for small and they can be neglected for small enough values of Q^. 
The numerical data for small indeed suggest a linear dependence. A linear 
extrapolation in Q^, using the three finest lattice spacings, gives us the values of the 
proton and neutron EDMs 


dp = 0.0340(62) 0 e-fm, (6.6) 

=-0.0318(54) d e-fm. (6.7) 

If we make the reasonable assumption that at this relatively large value of the 
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pseudoscalar mass, quenched and unquenched calculations give comparable results, 
we can try to estimate the values of the EDM at the physical point. To do so, we 
use as constraint the fact that the EDM in the continuum has to vanish in the 
chiral limit. In principle, we would like to use the yPT expressions in eqs. (2.4) and 
(2.5), but considering the large pseudoscalar mass used in our calculations, the yPT 
expressions are not reliable and we instead perform a simple linear £t in We 
then obtain the following estimates 

= 0.96(18) ■ 10“^ e e-im, ( 6 . 8 ) 

= -0.90(15) ■ 10"^ 9 e-fm, (6.9) 

where we have only included the errors from eqs. ( 6 . 6 ) and (6.7). These estimates are 
statistically consistent with the results of [12,13]. We stress that many systematic 
uncertainties are not taken into account in our calculation and that the main goal 
of this work is to describe the new methodology and perform a hrst continuum 
extrapolation. With all the caveats intrinsic in our calculation, we can extract an 
upper bound for 9. The experimental upper bound of the neutron EDM, |(ijv| < 
2.9 ■ 10-i3g. givgg 9 <3.2- 10-1°. 

Our lattice data also allows us to extract the nucleon Schiff moments for which we 
obtain 

=-1.16(33) ■ 10"^ 0 e ■ fm^ (6.10) 

= 1-07(28) ■ 10“^ 0 e-fm^. (6.11) 

In both SU{2) [18,6] and SU{3) [5,23] yPT, the nucleon Schiff moments are of 
isovector nature, in agreement with our lattice results. In fact, at leading order in 
the SU ( 2 ) chiral expansion, the Schiff moments are predicted 

■Sp = -SN = -jg^f|||^ = -(1.8±0,2).10-‘‘9e.fm=> , (6.12) 

where we used the value of given in eq. (2.9). 

We see that our lattice results are roughly 5 times larger than the leading-order yPT 
predictions which are in principle pion-mass independent (note that Qq ~ Tf^). How¬ 
ever, already at the physical pion mass, the nucleon Schiff moments obtain O{60%) 
next-to-leading-order corrections that scale as ~ [6,5]. Considering the large 

pseudoscalar mass used in our calculation, the discrepancy is not very worrisome. 

To conclude, we have shown that it is possible to obtain a non-perturbative deter¬ 
mination of the Schiff moments from the dependence of the CP-odd form factor, 
T 3 . With more precise lattice data at smaller pion masses, the extraction of the 
nucleon Schiff moments can be used for a direct determination of the EEC go. In 
this way, from the pion mass dependence of the EDMs (see eqs. (2.5,2.4)), it is then 
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possible to extract the other LECs, Ap^v (or equivalently d]\f^p), from lattice QCD 
alone without relying on eq. (2.9). 


7 Final remarks 


We have presented a first-of-its-kind continuum limit for the CP-mixing angle 
and the nucleon EDMs. At the same time, we have performed a first ab initio 
calculation of the nucleon Schiff moments in the continuum limit which can be 
used to extract the value of the CP-odd pion-nucleon coupling constant go. The key 
ingredient is the use of the gradient flow for the dehnition of the topological charge 
which, in this way, is free of renormalization ambiguities and allows a straightforward 
continuum extrapolation. The method we have proposed is general and can be used 
for any quantity computed in a 0 vacuum such as EDMs of light nuclei. To test 
this new method we have performed the calculation in the Yang-Mills theory at a 
relatively large value of the quark mass. We are currently extending this calculation 
to QCD with dynamical conhgurations and smaller quark masses. 

Previous calculations of nucleon EDMs [8,9], also applied a perturbative expansion 
in 9, but instead used a cooling procedure to dehne the topological charge. We 
do believe that such calculations give the right qualitative answer, but we stress 
that dehning the EDM in this way does not allow for a controlled continunm limit. 
Another method that has been proposed is to consider an imaginary 0-term in 
the action [14]. With an imaginary 6 the action becomes real and amenable to 
numerical Monte Carlo methods, bnt it requires the generation of a new gange 
ensemble for each value of 6. The range of 6 used for the generation of these ensembles 
is 0 ~ 0.5 — 2.5. If the calculation is performed with a Wilson-type action, the 6 
coefficient needs to be renormalized in order to restore the proper anomalous Ward 
identity. Additionally, the analytic continuation back to a real value of 6 has to be 
done with care in regions ontside the pertnrbative regime in 9. These complications 
are avoided completely with our proposal, because we directly compute the linear 
coefficient of 9 in the standard QCD backgronnd. 

As a last remark we recall that with a Wilson-type fermion action, EDMs are not 
guaranteed to vanish in the chiral limit, which only happens after the continunm 
limit has been performed. The same phenomenon takes place for the topological 
susceptibility [35]. It is only after performing the continuum limit that it is possible 
to constrain the 6'-induced EDM to vanish in the chiral limit. This stresses the 
importance of performing the continuum limit prior to any chiral limit and we 
believe that onr method is optimal in this respect. 

We consider this work as a hrst step in the determination of 9 and beyond-the- 
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Standard-Model matrix elements with the gradient flow. Other contributions to 
EDMs, for instance from fermionic operators such as quark chromo-electric dipole 
moments, can be determined using the gradient flow for fermions [28] and work in 
this direction is in progress. 
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A Appendix A 


In this appendix we discuss in some detail the spectral decomposition of the 3-point 
functions used to determine the nucleon form factors in a 6^ vacuum. For completeness 
we remind that in our calculations the initial momentum is pi = ((pi)o, Pi) and the 
final momentum is p2 = ((p 2 )o)P 2 )- The momentum transfer is q = p2 — Pi = 
(E(p2) - E(pi),P 2 - pi). 

If we retain only the leading exponential contribution, the spectral decomposition 
of the 3-point functions 

GT„a.(Pi.P2.2:o.!/o) = (A/'(x,i„)J,(y,s„)V( 0 ))^ , (A.l) 

x,y 

is given by 


^AfJ^v(Pl) P2; Xo, Vo) 


g-Ejv(pi)j/o g-Ejv(p2)(a:o-2/o) 
2EAr(pi) 2En{P2) 


(A.2) 


X j:{9|Af|iV''(p.2,s')>(iV''(P2.s')l4l'V'‘(Pi.s))('V'‘(Pi.s)|V|9). 


Following the parametrization in eq. (3.4) and using the completeness relation, for 
we obtain for small values of 6 
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g-i?jv(pi)j/o g—-E^jv(P2)(a^o—J/o) 

P 2 ; Xo, yo)ap = 2 En{p^) - 2Ejv(p2) —(A.3) 

X I [EN{p2ho - ilk{p2)k + Mn (1 + 2iea^^'’e'y^ 

X [EAr(pi)7o - i'ik{pi)k + Mm (1 + 2i6'aS^^6'75) . 

where al 3 are the Dirac indices. If we expand in powers of 6 the r.h.s. of eq. (A.l) 
we obtain 

G'wj^iv(Pn P 2 , xo, Vo) = Gnj^n{pu P 2 , xo, yo) + P 2 , xq, y ^), (A.4) 

where 

G'ivj^iv(Pi,P2,a;o,l/o) = ^A/'(x, xq) J^(y, l/o)A/'(0)) , (A.5) 

x,y 

is the 3-point function in the standard QCD background, and 

G'®j^,7v(Pi>P2,a^o,l/o) = a®^e*Piye*P2(^-y) ^7V'(x,a;o)J^(y,l/o)A/'(0)Q) , (A.6) 

x,y 

contains the insertion of the topological charge evaluated at non-vanishing flow-time 
■\/8t = O.Sro. 

Depending on the form factor we are interested in, we can select the appropriate 
Dirac indices with appropriate projectors that we indicate generically as 11 obtaining 

G'J(pi, P 2 ; xo, yo] n) = Tr [nG'^j^^(pi, pa; Xo, |/o)] , (A.7) 

i.e. 

Q-EN{Pl)yo g-Eiv(p 2 )(a:o-yo) 

GliPi,P2;xo,yo]n) = 2 Em{pi) - 2Ejv(p2) — 

xTr|n EAr(p 2 ) 7 o - i 7 fc(P 2 )fc + ATat (1-h 2i6'Q;j^^75) r^(Q^) 

X EAr(pi)7o - f7fc(pi)fc + Mm (1 + 2ieaM^'y5^ | . 

The spectral decomposition of the correlation functions in eqs. (A.5) and (A.6) 
traced with a generic projector 11 are easily obtained 

g-i?jv(pi)j/o g-Ejv(p2)(a:o-J/o) 

G'^(Pi,P2;a;o,|/o;n) = 2Em{pi) - 2Ejv(p2) — 

X Tr {n [E^(p 2 ) 7 o - t'yk{p 2 )k + Mm] r^™"(g2) [^iv(pi)7o - t7k{pi)k + Mm]} , 
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G'®(pi,P2;a;o,|/o;n) = 


g-i?jv(pi)j/0 g —-£^iv(P2)(a^0—J/o) 


Z*j^{Vi)Zn{V2) 


2-EAr(pi) 2i?7v(P2) 

X {Tr [n (2M^aff75) r7"(g') (^iv(pi)7o - *7fc(pi)fc + M^] 

n (^iv(P 2 ) 7 o - iiMk + Mn) r 7 "(g 2 ) ( 2 M^aff 75 ) 


(A.IO) 


Tr 


Tr 


n (-EAr(p2)7o — ilk{p2)k + Mm) r°‘^^(g^) (i?Ar(pi)7o — i')k{pi)k + Mm) } , 


where and r°^'^ are defined in eqs. (3.5) and (3.7). From this expression we 
already see that the 3-point fnnction with the insertion of the topological charge is 
not directly proportional to the CP-odd form factor but it contains additional 
contributions proportional to and the CP-even form factors. 

To extract the form factors traditionally one defines the following chain of ratios 

G'^(Pi,P2;a;o,|/o;n) 


^u(pi,P2;a;o,i/o;n) = 


C(p2,a;o) 

where C(p,a;o), the nucleon 2-point function, is defined as 

2 -£;jv(p)xo 


^(Pi,P2;a;o,|/o), 


C(p,a;o) = tr [P+GArjv(p,X q)] = \Zm{v)? ^ . (^7v(p) + Mm) + 

Fn{V) 


and 


^(Pi,P2;a;o,|/o) = 


C(P2, Xo)C(p2, |/o)(C(pi, Xq - Vo) 

C(pi,a;o)C(pi,|/o)C(P2,a;o - Vo) 


1/2 


(A.ll) 


(A.12) 


(A.13) 


For small 9 we have 


^J(Pi,P2;a;o,|/o;n) = P^(pi,p2;a;o,i/o;n) +i0Pj(pi,p2;a;o,|/o;n), (A.14) 


where 


D ( TT\ (pi ) P2 J 2^0) I/O) ff) 7^/ \ 

^p(Pi,P2;a;o,|/o;n) =--A:(pi,P2;a;o,|/o), 


and 


Ra{PuP2]Xo,yo-,U) = 


G^{Pi,P2-,xo,yo-,U) 

Cip2,Xo) 


■ K{pi,p2;xo,yo) . 


(A.15) 


(A.16) 


Performing the spectral decomposition and retaining only the fundamental state we 
obtain 



(A.17) 


72j(Pi, P 2 ; xo, vq] n) = A/'(pi, P 2 ) 

X Tr {n £^iv(p2)7o - hk{P 2 )k + Mat (1 + 2i9a^^^9j5 
^iv(Pi)7o - i'lk{Pi)k + Mn (1 + 2i6'aS^^6'75) } , 




where the normalization is given by 


A/'(pi,P2) = 


4-EAr(pi)-Ev(P2) 


_ -^Af(pi)-Ejv(p2) _ 

(EAr(pi) + M)(EAr(p2) + M) 


1/2 


(A.18) 


The ratio (A. 11) is dehned to remove the leading exponential contributions and have 
a plateau for 0 yo xq proportional to the form factors. 


We can now specialize the projector If, the external kinematics, and the current 
component in order to compute the form factors we want. For the 2 CP-even form 
factors we choose 


• n = P+, /i = 0, Pi = p, P2 = 0. 

In this case we indicate P(p) = E and P(p 2 ) = M 


Po(p, 0; xo, yo; P+) = A/'(p, 0) ■ 4MAr(PAr(p) -h Mn) 
where 

A/'(p,0) = 




4M2 


4M, 


N 


-|l/2 


2pAr(p)(Pv(p) + Mn) 


Putting everything together, we obtain 


Ro{p,0-,Xo,yo]P+) = 


En{p) + M]\f 
2Pv(p) 


1/2 




(A.19) 

(A.20) 

(A.21) 


To obtain eq. (A.21) and some of the eqs. below we have used the following 
kinematic relations 


|q|^ = |p|^ = Pl,-M^, 


(A.22) 


and 


= (Mat—Pat)^ —|q|^ ^ (f = (ATat—Pat)^—(P^—M^) = 2Mn{Mm—Em) < 0, 


This implies that 


Pat — Mtv — 


T 




(A.23) 


(A.24) 


2M]y 2M]\i 

Any of the relations in eq. (A.23) dehnes the = —q^ to be used when analyzing 
the form factors. 
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n = iP+757i, Pi = P, P2 = 0. 

After some algebra we obtain 


Ri{p, 0; xo, Vo; *P+757j) = A/'(p, 0) ■ AMg^eiju 


F.m+F^m 


(A.25) 


and using the expression for the normalization (A.20), we obtain 


Pi(p,0;a;o,i/o;*P+757i) 


1 


1/2 


2p7v(p)(pAf(p) + ATat) 


Qk^ijk 


F,{Q^) + F,m\ , 

(A.26) 


For the CP-odd form factor there are several choices for the Dirac projector. The 
analysis presented in this paper uses 


• n = iP+757i, h = 0, Pi = p, P2 = 0. 

If we compute the coefficient of i6, we obtain after some algebra 


RoiP, 0; xo, yo] *P+757i) = A/'(p, 0) ■ 


(A.27) 




Fim + 




2M 


N 


F2{Q^ 


2(Piv + M^)g,F3(Q')} 


Here the importance of a precise determination of becomes clear. The mixing of 
parity states induces spurious contributions to the correlation functions proportional 
to the CP-even form factors. These contributions need to be subtracted in order to 
determine the nucleon EDM. 
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